R Markdown

#Loading Packages
#Will most likely add more
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(bulletxtrctr)
## Warning in rgl.init(initValue, onlyNULL): RGL: unable to open X11 display
## Warning: 'rgl_init' failed, running with rgl.useNULL = TRUE
library(x3ptools)
library(randomForest)
## randomForest 4.6-12
## Type rfNews() to see new features/changes/bug fixes.
## 
## Attaching package: 'randomForest'
## The following object is masked from 'package:dplyr':
## 
##     combine
library(ggplot2)
## 
## Attaching package: 'ggplot2'
## The following object is masked from 'package:randomForest':
## 
##     margin
library(readr)
library(furrr)
## Loading required package: future
library(purrr)
library(stringr)
library(dichromat) 
library(future)
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## The following object is masked from 'package:randomForest':
## 
##     combine
## The following object is masked from 'package:dplyr':
## 
##     combine
library(tidyr)
options(future.globals.maxSize = 12*1024*1024*1024)
Group1 <- Group1 %>% 
  mutate(id = paste(Set, Barrel, Bullet, Land, sep = "-")) %>%
  select(id, Set, Barrel, Bullet, Land, dat, path) %>% 
  rename(x3p = dat)

Group2 <- Group2 %>% 
  mutate(id = paste(Set, Barrel, Bullet, Land, sep = "-")) %>%
  select(id, Set, Barrel, Bullet, Land, dat, path) %>% 
  rename(x3p = dat)

Group3 <- Group3 %>% 
  mutate(id = paste(Set, Barrel, Bullet, Land, sep = "-")) %>%
  select(id, Set, Barrel, Bullet, Land, dat, path) %>% 
  rename(x3p = dat)
Group1_Barrel_KA <- Group1 %>%
  filter(id %in% c("Group_1-KA- 1- 1", "Group_1-KA- 1- 2", "Group_1-KA- 1- 4", "Group_1-KA- 2- 1", "Group_1-KA- 2- 3", "Group_1-KA- 2- 5", "Group_1-KA- 2- 6", "Group_1-KA- 3- 1", "Group_1-KA- 3- 3", "Group_1-KA- 3- 4"))

Group1_Barrel_KB <- Group1 %>% 
  filter(id %in% c("Group_1-KB- 1- 2", "Group_1-KB- 1- 4", "Group_1-KB- 1- 6", "Group_1-KB- 2- 4", "Group_1-KB- 3- 1", "Group_1-KB- 3- 2"))


Group1_Barrel_KC <- Group1 %>% 
  filter(id %in% c("Group_1-KC- 1- 1", "Group_1-KC- 2- 3", "Group_1-KC- 2- 4", "Group_1-KC- 2- 5"))

Group1_Barrel_KD <- Group1 %>%
  filter(id %in% c("Group_1-KD- 1- 3", "Group_1-KD- 1- 4", "Group_1-KD- 2- 5", "Group_1-KD- 3- 1", "Group_1-KD- 3- 4"))

Group1_Barrel_KE <- Group1 %>% 
  filter(id %in% c("Group_1-KE- 2- 4"))

Group1_Barrel_Unknons <- Group1 %>%
  filter(Barrel == "Unknowns") 

#"Group_1-Unknowns-U37- 6", "Group_1-Unknowns-U42- 6", "Group_1-Unknowns-U77- 6", "Group_3-KJ- 3- 6"


Group1 <- Group1 %>%
  filter(!id %in% c("Group_1-KA- 1- 1", "Group_1-KA- 1- 2", "Group_1-KA- 1- 4", "Group_1-KA- 2- 1", "Group_1-KA- 2- 3", "Group_1-KA- 2- 5", "Group_1-KA- 2- 6", "Group_1-KA- 3- 1", "Group_1-KA- 3- 3", "Group_1-KA- 3- 4", "Group_1-KB- 1- 2", "Group_1-KB- 1- 4", "Group_1-KB- 1- 6", "Group_1-KB- 2- 4", "Group_1-KB- 3- 1", "Group_1-KB- 3- 2", "Group_1-KC- 1- 1", "Group_1-KC- 2- 3", "Group_1-KC- 2- 4", "Group_1-KC- 2- 5", "Group_1-KD- 1- 4", "Group_1-KD- 1- 3", "Group_1-KD- 2- 5", "Group_1-KD- 3- 1", "Group_1-KD- 3- 4", "Group_1-KE- 2- 4")) %>%
  filter(Barrel != "Unknowns")



plan(multicore) # use all the cores at once

#safe_crosscut <- safely(x3p_crosscut_optimize)

Group1 <- Group1 %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(150 ,NA)))

Group1_Barrel_KA <- Group1_Barrel_KA %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(350, NA)))

Group1_Barrel_KB <- Group1_Barrel_KB %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group1_Barrel_KC <- Group1_Barrel_KC %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group1_Barrel_KD <- Group1_Barrel_KD %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group1_Barrel_KE <- Group1_Barrel_KE %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group1_Barrel_Unknons <- Group1_Barrel_Unknons %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(250, NA)))


Group1 <- bind_rows(Group1, Group1_Barrel_KA, Group1_Barrel_KB, Group1_Barrel_KC, Group1_Barrel_KD, Group1_Barrel_KE, Group1_Barrel_Unknons)
Group2_Barrel_KC <- Group2 %>%
  filter(id %in% c("Group_2-KC- 1- 3", "Group_2-KC- 1- 4", "Group_2-KC- 2- 2", "Group_2-KC- 2- 5", "Group_2-KC- 2- 6", "Group_2-KC- 3- 2", "Group_2-KC- 3- 6"))

Group2_Barrel_KD <- Group2 %>% 
  filter(id %in% c("Group_2-KD- 1- 1", "Group_2-KD- 1- 2", "Group_2-KD- 1- 3", "Group_2-KD- 1- 4", "Group_2-KD- 3- 2", "Group_2-KD- 3- 3", "Group_2-KD- 3- 4", "Group_2-KD- 3- 5"))

Group2_Barrel_KE <- Group2 %>% 
  filter(id %in% c("Group_2-KE- 1- 3", "Group_2-KE- 1- 5", "Group_2-KE- 2- 2", "Group_2-KE- 2- 3", "Group_2-KE- 3- 2", "Group_2-KE- 3- 5"))

Group2_Barrel_KF <- Group2 %>% 
  filter(id %in% c("Group_2-KF- 1- 2", "Group_2-KF- 2- 1", "Group_2-KF- 2- 3", "Group_2-KF- 2- 5", "Group_2-KF- 2- 6", "Group_2-KF- 3- 3"))

Group2_Barrel_KG <- Group2 %>%
  filter(id %in% c("Group_2-KG- 1- 1", "Group_2-KG- 1- 5", "Group_2-KG- 1- 6", "Group_2-KG- 2- 3", "Group_2-KG- 2- 4", "Group_2-KG- 3- 3"))

Group2_Barrel_Unknowns <- Group2 %>%
  filter(id %in% c("Group_2-Unknowns-U23- 1", "Group_2-Unknowns-U23- 3", "Group_2-Unknowns-U23- 6", "Group_2-Unknowns-U41- 1", "Group_2-Unknowns-U41- 2", "Group_2-Unknowns-U61- 4", "Group_2-Unknowns-U63- 2", "Group_2-Unknowns-U63- 4", "Group_2-Unknowns-U66- 2", "Group_2-Unknowns-U66- 5", "Group_2-Unknowns-U73- 2", "Group_2-Unknowns-U73- 5", "Group_2-Unknowns-U98- 4", "Group_2-Unknowns-U98- 6"))

Group2 <- Group2 %>%
  filter(!id %in% c("Group_2-KC- 1- 3", "Group_2-KC- 1- 4", "Group_2-KC- 2- 2", "Group_2-KC- 2- 5", "Group_2-KC- 2- 6", "Group_2-KC- 3- 2", "Group_2-KC- 3- 6", "Group_2-KD- 1- 1", "Group_2-KD- 1- 2", "Group_2-KD- 1- 3", "Group_2-KD- 1- 4", "Group_2-KD- 3- 2", "Group_2-KD- 3- 3", "Group_2-KD- 3- 4", "Group_2-KD- 3- 5", "Group_2-KE- 1- 3", "Group_2-KE- 1- 5", "Group_2-KE- 2- 2", "Group_2-KE- 2- 3", "Group_2-KE- 3- 2", "Group_2-KE- 3- 5", "Group_2-KF- 1- 2", "Group_2-KF- 2- 1", "Group_2-KF- 2- 3", "Group_2-KF- 2- 5", "Group_2-KF- 2- 6", "Group_2-KF- 3- 3", "Group_2-KG- 1- 1", "Group_2-KG- 1- 5", "Group_2-KG- 1- 6", "Group_2-KG- 2- 3", "Group_2-KG- 2- 4", "Group_2-KG- 3- 3", "Group_2-Unknowns-U23- 1", "Group_2-Unknowns-U23- 3", "Group_2-Unknowns-U23- 6", "Group_2-Unknowns-U41- 1", "Group_2-Unknowns-U41- 2", "Group_2-Unknowns-U61- 4", "Group_2-Unknowns-U63- 2", "Group_2-Unknowns-U63- 4", "Group_2-Unknowns-U66- 2", "Group_2-Unknowns-U66- 5", "Group_2-Unknowns-U73- 2", "Group_2-Unknowns-U73- 5", "Group_2-Unknowns-U98- 4", "Group_2-Unknowns-U98- 6"))



Group2 <- Group2 %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize))

Group2_Barrel_KC <- Group2_Barrel_KC %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group2_Barrel_KD <- Group2_Barrel_KD %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group2_Barrel_KE <- Group2_Barrel_KE %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group2_Barrel_KF <- Group2_Barrel_KF %>% 
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(300, NA)))

Group2_Barrel_KG <- Group2_Barrel_KG %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(350, NA)))

Group2_Barrel_Unknowns <- Group2_Barrel_Unknowns %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimits = c(350, NA)))

Group2 <- bind_rows(Group2, Group2_Barrel_KC, Group2_Barrel_KD, Group2_Barrel_KE, Group2_Barrel_KF, Group2_Barrel_KG, Group2_Barrel_Unknowns)
Group3 <- Group3 %>%
  mutate(CrossSection = future_map_dbl(x3p, x3p_crosscut_optimize, ylimit = c(250, NA)))
Group1 <- Group1 %>% 
  mutate(CrossCut = future_map2(.x = x3p, .y = CrossSection, .f = x3p_crosscut))

Houston_CrossCuts_1 <- select(Group1, -path, -x3p) %>% 
      tidyr::unnest(CrossCut)
Group2 <- Group2 %>% 
  mutate(CrossCut = future_map2(.x = x3p, .y = CrossSection, .f = x3p_crosscut))

Houston_CrossCuts_2 <- select(Group2, -path, -x3p) %>%
  tidyr::unnest(CrossCut)
Group3 <- Group3 %>% 
  mutate(CrossCut = future_map2(.x = x3p, .y = CrossSection, .f = x3p_crosscut))

Houston_CrossCuts_3 <- select(Group3, -path, -x3p) %>%
  tidyr::unnest(CrossCut)
Houston_CrossCuts_1 %>%
  filter(Barrel != "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_grid(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group1 Known Crosscut values")

Houston_CrossCuts_1 %>%
  filter(Barrel == "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_wrap(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group1 Unknown Crosscut values")

Houston_CrossCuts_2 %>%
  filter(Barrel != "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_grid(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group2 Known Crosscut values")

Houston_CrossCuts_2 %>%
  filter(Barrel == "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_wrap(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group2 Unknown Crosscut values")

Houston_CrossCuts_3 %>%
  filter(Barrel != "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_grid(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group3 Known Crosscut values")

Houston_CrossCuts_3 %>%
  filter(Barrel == "Unknowns") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_wrap(Set + paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group3 Unknown Crosscut values")

Houston_CrossCuts_3 %>%
  filter(Barrel == "KJ") %>% 
    ggplot(data = ., aes(x = x, y = value, color = Land)) + 
      geom_line() + 
        facet_wrap(~Bullet+Land, ncol = 6)+
          theme_bw()

#checking for possible better regions for crosscut scans
head(Group1, 138)
## # A tibble: 138 x 9
##    id     Set    Barrel Bullet Land  x3p   path       CrossSection CrossCut
##    <chr>  <chr>  <chr>  <chr>  <chr> <lis> <chr>             <dbl> <list>  
##  1 Group… Group… KA     " 1"   " 3"  <x3p> /media/Su…          150 <df[,3]…
##  2 Group… Group… KA     " 1"   " 5"  <x3p> /media/Su…          150 <df[,3]…
##  3 Group… Group… KA     " 1"   " 6"  <x3p> /media/Su…          150 <df[,3]…
##  4 Group… Group… KA     " 2"   " 2"  <x3p> /media/Su…          150 <df[,3]…
##  5 Group… Group… KA     " 2"   " 4"  <x3p> /media/Su…          150 <df[,3]…
##  6 Group… Group… KA     " 3"   " 2"  <x3p> /media/Su…          150 <df[,3]…
##  7 Group… Group… KA     " 3"   " 5"  <x3p> /media/Su…          150 <df[,3]…
##  8 Group… Group… KA     " 3"   " 6"  <x3p> /media/Su…          150 <df[,3]…
##  9 Group… Group… KB     " 1"   " 1"  <x3p> /media/Su…          150 <df[,3]…
## 10 Group… Group… KB     " 1"   " 3"  <x3p> /media/Su…          150 <df[,3]…
## # … with 128 more rows
head(Group2, 138)
## # A tibble: 138 x 9
##    id     Set    Barrel Bullet Land  x3p   path       CrossSection CrossCut
##    <chr>  <chr>  <chr>  <chr>  <chr> <lis> <chr>             <dbl> <list>  
##  1 Group… Group… KC     " 1"   " 1"  <x3p> /media/Su…          100 <df[,3]…
##  2 Group… Group… KC     " 1"   " 2"  <x3p> /media/Su…           75 <df[,3]…
##  3 Group… Group… KC     " 1"   " 5"  <x3p> /media/Su…          150 <df[,3]…
##  4 Group… Group… KC     " 1"   " 6"  <x3p> /media/Su…          175 <df[,3]…
##  5 Group… Group… KC     " 2"   " 1"  <x3p> /media/Su…          125 <df[,3]…
##  6 Group… Group… KC     " 2"   " 3"  <x3p> /media/Su…          100 <df[,3]…
##  7 Group… Group… KC     " 2"   " 4"  <x3p> /media/Su…           75 <df[,3]…
##  8 Group… Group… KC     " 3"   " 1"  <x3p> /media/Su…          150 <df[,3]…
##  9 Group… Group… KC     " 3"   " 3"  <x3p> /media/Su…          150 <df[,3]…
## 10 Group… Group… KC     " 3"   " 4"  <x3p> /media/Su…          125 <df[,3]…
## # … with 128 more rows
head(Group3, 138)
## # A tibble: 138 x 9
##    id     Set    Barrel Bullet Land  x3p   path       CrossSection CrossCut
##    <chr>  <chr>  <chr>  <chr>  <chr> <lis> <chr>             <dbl> <list>  
##  1 Group… Group… KF     " 1"   " 1"  <x3p> /media/Su…          250 <df[,3]…
##  2 Group… Group… KF     " 1"   " 2"  <x3p> /media/Su…          250 <df[,3]…
##  3 Group… Group… KF     " 1"   " 3"  <x3p> /media/Su…          325 <df[,3]…
##  4 Group… Group… KF     " 1"   " 4"  <x3p> /media/Su…          275 <df[,3]…
##  5 Group… Group… KF     " 1"   " 5"  <x3p> /media/Su…          250 <df[,3]…
##  6 Group… Group… KF     " 1"   " 6"  <x3p> /media/Su…          250 <df[,3]…
##  7 Group… Group… KF     " 2"   " 1"  <x3p> /media/Su…          250 <df[,3]…
##  8 Group… Group… KF     " 2"   " 2"  <x3p> /media/Su…          300 <df[,3]…
##  9 Group… Group… KF     " 2"   " 3"  <x3p> /media/Su…          250 <df[,3]…
## 10 Group… Group… KF     " 2"   " 4"  <x3p> /media/Su…          250 <df[,3]…
## # … with 128 more rows
# Grooves 

saved_grooves_location_Houston_1 <- "Group1data.rda"
if (file.exists(saved_grooves_location_Houston_1)) {
  Group1$Grooves <- readRDS(saved_grooves_location_Houston_1)
} else {
  Group1 <- Group1 %>% 
    mutate(Grooves = CrossCut %>% 
             future_map(.f = cc_locate_grooves, 
                        method = "rollapply", smoothfactor = 15, return_plot = T))  # use plot so that the shiny app works...
}

grooves_Group1 <- Group1 %>% tidyr::unnest(Grooves)
# Grooves 

saved_grooves_location_Houston_2 <- "Group2data.rda"
if (file.exists(saved_grooves_location_Houston_2)) {
  Group2$Grooves <- readRDS(saved_grooves_location_Houston_2)
} else {
  Group2 <- Group2 %>% 
    mutate(Grooves = CrossCut %>% 
             future_map(.f = cc_locate_grooves, 
                        method = "rollapply", smoothfactor = 15, return_plot = T))  # use plot so that the shiny app works...
}

grooves_Group2 <- Group2 %>% tidyr::unnest(Grooves)
# Grooves

saved_grooves_location_Houston_3 <- "Group3data.rda"
if (file.exists(saved_grooves_location_Houston_3)) {
  Group3$Grooves <- readRDS(saved_grooves_location_Houston_3)
} else {
  Group3 <- Group3 %>% 
    mutate(Grooves = CrossCut %>% 
             future_map(.f = cc_locate_grooves, 
                        method = "rollapply", smoothfactor = 15, return_plot = T))  # use plot so that the shiny app works...
}

grooves_Group3 <- Group3 %>% tidyr::unnest(Grooves)

Shiny app

library(shiny)
if (file.exists(saved_grooves_location_Houston_1)) {
  Group1$Grooves <- readRDS(saved_grooves_location_Houston_1)
}
if (interactive()) { # only run when you're manually running chunks... don't run when the whole document is compiled.
  shinyApp(
    ui = fluidPage(
      selectInput("k", "Investigate kth plot:",
                  selected = 1,
                  choices = (1:length(Group1$Grooves)) %>% set_names(Group1$id)
      ),
      textOutput("groovelocations"),
      actionButton("confirm", "Confirm"),
      actionButton("save", "Save"),
      plotOutput("groovePlot", click = "plot_click"),
      verbatimTextOutput("info")
    ),
    
    server = function(input, output, session) {
      output$groovePlot <- renderPlot({
        k <- as.numeric(input$k)
        p <- Group1$Grooves[[k]]$plot
        
        p
      })
      output$groovelocations <- renderText({
        paste(
          "Left Groove: ", Group1$Grooves[[as.numeric(input$k)]]$groove[1],
          " Right Groove: ", Group1$Grooves[[as.numeric(input$k)]]$groove[2]
        )
      })
      observeEvent(input$confirm, {
        cat(paste(Group1$id[as.numeric(input$k)], "\n"))
        updateSelectInput(session, "k", "Investigate kth plot:",
                          selected = as.numeric(input$k) + 1,
                          choices = (1:length(Group1$Grooves)) %>% set_names(Group1$id)
        )
      })
      observeEvent(input$save, {
        saveRDS(Group1$Grooves, file = saved_grooves_location_Houston_1)
        message("groove data saved\n")
      })
      
      observeEvent(input$plot_click, {
        k <- as.numeric(input$k)
        xloc <- input$plot_click$x
        
        gr <- Group1$Grooves[[k]]$groove
        if (abs(gr[1] - xloc) < abs(gr[2] - xloc)) {
          Group1$Grooves[[k]]$groove[1] <<- xloc
        } else {
          Group1$Grooves[[k]]$groove[2] <<- xloc
        }
        output$groovePlot <- renderPlot({
          k <- as.numeric(input$k)
          p <- Group1$Grooves[[k]]$plot +
            geom_vline(xintercept = Group1$Grooves[[k]]$groove[1], colour = "green") +
            geom_vline(xintercept = Group1$Grooves[[k]]$groove[2], colour = "green")
          
          p
        })
      })
      output$info <- renderText({
        paste0("x=", input$plot_click$x, "\ny=", input$plot_click$y)
      })
    },
    options = list(height = 500)
  )
  saveRDS(Group1$Grooves, file = saved_grooves_location_Houston_1)
} else {
  if (!file.exists(saved_grooves_location_Houston_1)) {
    message("run script in interactive mode to fix grooves")
  } else {
    Group1$Grooves <- readRDS(saved_grooves_location_Houston_1)
  }
}
library(shiny)
if (file.exists(saved_grooves_location_Houston_2)) {
  Group2$Grooves <- readRDS(saved_grooves_location_Houston_2)
}
if (interactive()) { # only run when you're manually running chunks... don't run when the whole document is compiled.
  shinyApp(
    ui = fluidPage(
      selectInput("k", "Investigate kth plot:",
                  selected = 1,
                  choices = (1:length(Group2$Grooves)) %>% set_names(Group2$id)
      ),
      textOutput("groovelocations"),
      actionButton("confirm", "Confirm"),
      actionButton("save", "Save"),
      plotOutput("groovePlot", click = "plot_click"),
      verbatimTextOutput("info")
    ),
    
    server = function(input, output, session) {
      output$groovePlot <- renderPlot({
        k <- as.numeric(input$k)
        p <- Group2$Grooves[[k]]$plot
        
        p
      })
      output$groovelocations <- renderText({
        paste(
          "Left Groove: ", Group2$Grooves[[as.numeric(input$k)]]$groove[1],
          " Right Groove: ", Group2$Grooves[[as.numeric(input$k)]]$groove[2]
        )
      })
      observeEvent(input$confirm, {
        cat(paste(Group2$id[as.numeric(input$k)], "\n"))
        updateSelectInput(session, "k", "Investigate kth plot:",
                          selected = as.numeric(input$k) + 1,
                          choices = (1:length(Group2$Grooves)) %>% set_names(Group2$id)
        )
      })
      observeEvent(input$save, {
        saveRDS(Group2$Grooves, file = saved_grooves_location_Houston_2)
        message("groove data saved\n")
      })
      
      observeEvent(input$plot_click, {
        k <- as.numeric(input$k)
        xloc <- input$plot_click$x
        
        gr <- Group2$Grooves[[k]]$groove
        if (abs(gr[1] - xloc) < abs(gr[2] - xloc)) {
          Group2$Grooves[[k]]$groove[1] <<- xloc
        } else {
          Group2$Grooves[[k]]$groove[2] <<- xloc
        }
        output$groovePlot <- renderPlot({
          k <- as.numeric(input$k)
          p <- Group2$Grooves[[k]]$plot +
            geom_vline(xintercept = Group2$Grooves[[k]]$groove[1], colour = "green") +
            geom_vline(xintercept = Group2$Grooves[[k]]$groove[2], colour = "green")
          
          p
        })
      })
      output$info <- renderText({
        paste0("x=", input$plot_click$x, "\ny=", input$plot_click$y)
      })
    },
    options = list(height = 500)
  )
  saveRDS(Group2$Grooves, file = saved_grooves_location_Houston_2)
} else {
  if (!file.exists(saved_grooves_location_Houston_2)) {
    message("run script in interactive mode to fix grooves")
  } else {
    Group2$Grooves <- readRDS(saved_grooves_location_Houston_2)
  }
}
library(shiny)
if (file.exists(saved_grooves_location_Houston_3)) {
  Group3$Grooves <- readRDS(saved_grooves_location_Houston_3)
}
if (interactive()) { # only run when you're manually running chunks... don't run when the whole document is compiled.
  shinyApp(
    ui = fluidPage(
      selectInput("k", "Investigate kth plot:",
                  selected = 1,
                  choices = (1:length(Group3$Grooves)) %>% set_names(Group3$id)
      ),
      textOutput("groovelocations"),
      actionButton("confirm", "Confirm"),
      actionButton("save", "Save"),
      plotOutput("groovePlot", click = "plot_click"),
      verbatimTextOutput("info")
    ),
    
    server = function(input, output, session) {
      output$groovePlot <- renderPlot({
        k <- as.numeric(input$k)
        p <- Group3$Grooves[[k]]$plot
        
        p
      })
      output$groovelocations <- renderText({
        paste(
          "Left Groove: ", Group3$Grooves[[as.numeric(input$k)]]$groove[1],
          " Right Groove: ", Group3$Grooves[[as.numeric(input$k)]]$groove[2]
        )
      })
      observeEvent(input$confirm, {
        cat(paste(Group3$id[as.numeric(input$k)], "\n"))
        updateSelectInput(session, "k", "Investigate kth plot:",
                          selected = as.numeric(input$k) + 1,
                          choices = (1:length(Group3$Grooves)) %>% set_names(Group3$id)
        )
      })
      observeEvent(input$save, {
        saveRDS(Group3$Grooves, file = saved_grooves_location_Houston_3)
        message("groove data saved\n")
      })
      
      observeEvent(input$plot_click, {
        k <- as.numeric(input$k)
        xloc <- input$plot_click$x
        
        gr <- Group3$Grooves[[k]]$groove
        if (abs(gr[1] - xloc) < abs(gr[2] - xloc)) {
          Group3$Grooves[[k]]$groove[1] <<- xloc
        } else {
          Group3$Grooves[[k]]$groove[2] <<- xloc
        }
        output$groovePlot <- renderPlot({
          k <- as.numeric(input$k)
          p <- Group3$Grooves[[k]]$plot +
            geom_vline(xintercept = Group3$Grooves[[k]]$groove[1], colour = "green") +
            geom_vline(xintercept = Group3$Grooves[[k]]$groove[2], colour = "green")
          
          p
        })
      })
      output$info <- renderText({
        paste0("x=", input$plot_click$x, "\ny=", input$plot_click$y)
      })
    },
    options = list(height = 500)
  )
  saveRDS(Group3$Grooves, file = saved_grooves_location_Houston_3)
} else {
  if (!file.exists(saved_grooves_location_Houston_3)) {
    message("run script in interactive mode to fix grooves")
  } else {
    Group3$Grooves <- readRDS(saved_grooves_location_Houston_3)
  }
}
Group1 <- Group1 %>% 
 mutate(Signatures = future_map2(.x = CrossCut, .y = Grooves, .f = cc_get_signature, span = 0.75, span2 = .03))

 Signatures_Group1 <- Group1 %>% 
  select(id, Set, Barrel, Bullet, Land, Signatures) %>% 
   tidyr::unnest()

Signatures_Group1 %>%
  filter(Barrel != "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_grid(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group1 Signatures Known")
## Warning: Removed 4983 rows containing missing values (geom_path).

Signatures_Group1 %>%
  filter(Barrel == "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_wrap(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group1 Signatures Unknown")
## Warning: Removed 6386 rows containing missing values (geom_path).

Group2 <- Group2 %>% 
 mutate(Signatures = future_map2(.x = CrossCut, .y = Grooves, .f = cc_get_signature, span = 0.75, span2 = .03))

 Signatures_Group2 <- Group2 %>% 
  select(id, Set, Barrel, Bullet, Land, Signatures) %>% 
   tidyr::unnest()

Signatures_Group2 %>%
  filter(Barrel != "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_grid(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group2 Signatures Known")
## Warning: Removed 5359 rows containing missing values (geom_path).

Signatures_Group2 %>%
  filter(Barrel == "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_wrap(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group2 Signatures Unknown")
## Warning: Removed 5770 rows containing missing values (geom_path).

Group3 <- Group3 %>% 
 mutate(Signatures = future_map2(.x = CrossCut, .y = Grooves, .f = cc_get_signature, span = 0.75, span2 = .03))

 Signatures_Group3 <- Group3 %>% 
  select(id, Set, Barrel, Bullet, Land, Signatures) %>% 
   tidyr::unnest()

Signatures_Group3 %>%
  filter(Barrel != "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_grid(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group3 Signatures Known")
## Warning: Removed 4566 rows containing missing values (geom_path).

Signatures_Group3 %>%
  filter(Barrel == "Unknowns") %>%
    ggplot(data = ., aes(x = x, y = sig, color = Land)) + 
      geom_line()+
        facet_wrap(paste("Bullet", Bullet) ~ sprintf("Barrel %02s", Barrel))+
          theme_bw()+
            ggtitle("Group3 Signatures Unknown")
## Warning: Removed 6077 rows containing missing values (geom_path).

comparisons_1 <- crossing(Bullet1 = Group1$id, Bullet2 = Group1$id) %>%
  left_join(nest(Group1, -id) %>% magrittr::set_names(c("Bullet1", "Bullet1_data"))) %>%
  left_join(nest(Group1, -id) %>% magrittr::set_names(c("Bullet2", "Bullet2_data"))) %>%
  mutate(Set1 = str_extract(Bullet1, "Group_\\d{1}"),
         Set2 = str_extract(Bullet2, "Group_\\d{1}")) %>%
  filter(Set1 == Set2) %>% # Get rid of cross-set comparisons for now...
  select(-matches("Set"))
## Joining, by = "Bullet1"
## Joining, by = "Bullet2"
comparisons_2 <- crossing(Bullet1 = Group2$id, Bullet2 = Group2$id) %>%
  left_join(nest(Group2, -id) %>% magrittr::set_names(c("Bullet1", "Bullet1_data"))) %>%
  left_join(nest(Group2, -id) %>% magrittr::set_names(c("Bullet2", "Bullet2_data"))) %>%
  mutate(Set1 = str_extract(Bullet1, "Group_\\d{1}"),
         Set2 = str_extract(Bullet2, "Group_\\d{1}")) %>%
  filter(Set1 == Set2) %>% # Get rid of cross-set comparisons for now...
  select(-matches("Set"))
## Joining, by = "Bullet1"
## Joining, by = "Bullet2"
comparisons_3 <- crossing(Bullet1 = Group3$id, Bullet2 = Group3$id) %>%
  left_join(nest(Group3, -id) %>% magrittr::set_names(c("Bullet1", "Bullet1_data"))) %>%
  left_join(nest(Group3, -id) %>% magrittr::set_names(c("Bullet2", "Bullet2_data"))) %>%
  mutate(Set1 = str_extract(Bullet1, "Group_\\d{1}"),
         Set2 = str_extract(Bullet2, "Group_\\d{1}")) %>%
  filter(Set1 == Set2) %>% # Get rid of cross-set comparisons for now...
  select(-matches("Set"))
## Joining, by = "Bullet1"
## Joining, by = "Bullet2"
#plan(multicore(workers = availableCores(constraints = 8)))

plan(multicore)

get_sig <- function(data) {
  map(data$Signatures, "sig")
}
comparisons_1 <- comparisons_1 %>%
  mutate(sig1 = future_map(Bullet1_data, get_sig), sig2 = future_map(Bullet2_data, get_sig))

comparisons_1 <- comparisons_1 %>%
  mutate(Aligned = future_map2(sig1, sig2, ~sig_align(unlist(.x), unlist(.y)))) # Getting Aligned signatures
# Get striae

comparisons_1 <- comparisons_1 %>%
  mutate(Striae = future_map(Aligned, sig_cms_max)) # Obtaining Striae
saveRDS(select(comparisons_1, -Bullet1_data, -Bullet2_data), file = "Group1_Comparisons.rda")


comparisons_1 <- comparisons_1 %>% 
  select(-Bullet1_data, -Bullet2_data)

comparisons_1 <- comparisons_1 %>% 
  mutate(features = future_map2(.x = Aligned, .y = Striae, .f = extract_features_all, resolution = 1.5625))#ObtainingFeatures

comparisons_1 <- comparisons_1 %>% 
  mutate(Legacy_Features = future_map(Striae, extract_features_all_legacy, resolution = 1.5625)) # Obtaining feature leacy

comparisons_legacy_1 <- comparisons_1 %>% 
  select(-features) %>%
  tidyr::unnest(Legacy_Features) # Extracting feature legacy

comparisons_1 <- comparisons_1 %>%
  select(-Legacy_Features) %>%
  tidyr::unnest(features)
#plan(multicore(workers = availableCores(constraints = 8)))

plan(multicore)

get_sig <- function(data) {
  map(data$Signatures, "sig")
}
comparisons_2 <- comparisons_2 %>%
  mutate(sig1 = future_map(Bullet1_data, get_sig), sig2 = future_map(Bullet2_data, get_sig))

comparisons_2 <- comparisons_2 %>%
  mutate(Aligned = future_map2(sig1, sig2, ~sig_align(unlist(.x), unlist(.y)))) # Getting Aligned signatures
# Get striae

comparisons_2 <- comparisons_2 %>%
  mutate(Striae = future_map(Aligned, sig_cms_max)) # Obtaining Striae
saveRDS(select(comparisons_2, -Bullet1_data, -Bullet2_data), file = "Group2_Comparisons.rda")


comparisons_2 <- comparisons_2 %>% 
  select(-Bullet1_data, -Bullet2_data)

comparisons_2 <- comparisons_2 %>% 
  mutate(features = future_map2(.x = Aligned, .y = Striae, .f = extract_features_all, resolution = 1.5625))#ObtainingFeatures

comparisons_2 <- comparisons_2 %>% 
  mutate(Legacy_Features = future_map(Striae, extract_features_all_legacy, resolution = 1.5625)) # Obtaining feature leacy

comparisons_legacy_2 <- comparisons_2 %>% 
  select(-features) %>%
  tidyr::unnest(Legacy_Features) # Extracting feature legacy

comparisons_2 <- comparisons_2 %>%
  select(-Legacy_Features) %>%
  tidyr::unnest(features)
plan(multicore)

get_sig <- function(data) {
  map(data$Signatures, "sig")
}
comparisons_3 <- comparisons_3 %>%
  mutate(sig1 = future_map(Bullet1_data, get_sig), sig2 = future_map(Bullet2_data, get_sig))

comparisons_3 <- comparisons_3 %>%
  mutate(Aligned = future_map2(sig1, sig2, ~sig_align(unlist(.x), unlist(.y)))) # Getting Aligned signatures
# Get striae

comparisons_3 <- comparisons_3 %>%
  mutate(Striae = future_map(Aligned, sig_cms_max)) # Obtaining Striae
saveRDS(select(comparisons_3, -Bullet1_data, -Bullet2_data), file = "Group3_Comparisons.rda")



comparisons_3 <- comparisons_3 %>% 
  select(-Bullet1_data, -Bullet2_data)

comparisons_3 <- comparisons_3 %>% 
  mutate(features = future_map2(.x = Aligned, .y = Striae, .f = extract_features_all, resolution = 1.5625))#ObtainingFeatures

comparisons_3 <- comparisons_3 %>% 
  mutate(Legacy_Features = future_map(Striae, extract_features_all_legacy, resolution = 1.5625)) # Obtaining feature leacy

comparisons_legacy_3 <- comparisons_3 %>% 
  select(-features) %>%
  tidyr::unnest(Legacy_Features) # Extracting feature legacy

comparisons_3 <- comparisons_3 %>%
  #select(-Legacy_Features) %>%
  tidyr::unnest(features)

head(comparisons_1)
## # A tibble: 6 x 28
##   Bullet1 Bullet2 sig1  sig2  Aligned Striae   ccf   cms  cms2 cms2_per_mm
##   <chr>   <chr>   <lis> <lis> <list>  <list> <dbl> <dbl> <dbl>       <dbl>
## 1 Group_… Group_… <lis… <lis… <named… <name… 1        27    14       4.74 
## 2 Group_… Group_… <lis… <lis… <named… <name… 0.566     1     2       0.679
## 3 Group_… Group_… <lis… <lis… <named… <name… 0.474     3     3       1.03 
## 4 Group_… Group_… <lis… <lis… <named… <name… 0.552     2     2       0.677
## 5 Group_… Group_… <lis… <lis… <named… <name… 0.400     0     0       0    
## 6 Group_… Group_… <lis… <lis… <named… <name… 0.385     5     4       1.35 
## # … with 18 more variables: cms_per_mm <dbl>, D <dbl>, lag <dbl>,
## #   lag_mm <dbl>, left_cms <dbl>, length <dbl>, length_mm <dbl>,
## #   matches <dbl>, matches_per_mm <dbl>, mismatches <dbl>,
## #   mismatches_per_mm <dbl>, non_cms <dbl>, non_cms_per_mm <dbl>,
## #   overlap <dbl>, right_cms <dbl>, rough_cor <dbl>, sd_D <dbl>,
## #   sum_peaks <dbl>
comparisons_1 <- comparisons_1 %>%
  select(-sig1, -sig2, -Aligned, -Striae, - right_cms, -left_cms)

comparisons_1 <- comparisons_1 %>%
  mutate(Bullet1 = gsub(" ", "", fixed = TRUE, Bullet1), Bullet2 = gsub(" ", "", fixed = TRUE, Bullet2))

comparisons_1 <- comparisons_1 %>% 
  mutate(Set = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\1", Bullet2))
comparisons_1 <- comparisons_1 %>% 
  mutate(BarrelA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet2))
comparisons_1 <- comparisons_1 %>% 
  mutate(BarrelB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet1))
comparisons_1 <- comparisons_1 %>% 
  mutate(BulletA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet2))
comparisons_1 <- comparisons_1 %>% 
  mutate(BulletB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet1))
comparisons_1 <- comparisons_1 %>% 
  mutate(LandA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet2))
comparisons_1 <- comparisons_1 %>% 
  mutate(LandB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet1))
comparisons_2 <- comparisons_2 %>%
  select(-sig1, -sig2, -Aligned, -Striae, - right_cms, -left_cms)

comparisons_2 <- comparisons_2 %>%
  mutate(Bullet1 = gsub(" ", "", fixed = TRUE, Bullet1), Bullet2 = gsub(" ", "", fixed = TRUE, Bullet2))

comparisons_2 <- comparisons_2 %>% 
  mutate(Set = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\1", Bullet2))
comparisons_2 <- comparisons_2 %>% 
  mutate(BarrelA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet2))
comparisons_2 <- comparisons_2 %>% 
  mutate(BarrelB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet1))
comparisons_2 <- comparisons_2 %>% 
  mutate(BulletA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet2))
comparisons_2 <- comparisons_2 %>% 
  mutate(BulletB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet1))
comparisons_2 <- comparisons_2 %>% 
  mutate(LandA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet2))
comparisons_2 <- comparisons_2 %>% 
  mutate(LandB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet1))
comparisons_3 <- comparisons_3 %>%
  select(-sig1, -sig2, -Aligned, -Striae, - right_cms, -left_cms)

comparisons_3 <- comparisons_3 %>%
  mutate(Bullet1 = gsub(" ", "", fixed = TRUE, Bullet1), Bullet2 = gsub(" ", "", fixed = TRUE, Bullet2))

comparisons_3 <- comparisons_3 %>% 
  mutate(Set = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\1", Bullet2))
comparisons_3 <- comparisons_3 %>% 
  mutate(BarrelA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet2))
comparisons_3 <- comparisons_3 %>% 
  mutate(BarrelB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\2", Bullet1))
comparisons_3 <- comparisons_3 %>% 
  mutate(BulletA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet2))
comparisons_3 <- comparisons_3 %>% 
  mutate(BulletB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\3", Bullet1))
comparisons_3 <- comparisons_3 %>% 
  mutate(LandA = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet2))
comparisons_3 <- comparisons_3 %>% 
  mutate(LandB = gsub("(Group_[0-9]{1})-([A-Z]{1,2}|Unknowns)-(U[0-9]{1,2}|[0-9A-Z])-([1-6])", "\\4", Bullet1))

head(comparisons_2)
## # A tibble: 6 x 29
##   Bullet1 Bullet2   ccf   cms  cms2 cms2_per_mm cms_per_mm      D   lag
##   <chr>   <chr>   <dbl> <dbl> <dbl>       <dbl>      <dbl>  <dbl> <dbl>
## 1 Group_… Group_… 1        27    13       4.71       9.78  0          0
## 2 Group_… Group_… 0.536     1     0       0          0.362 0.0656    52
## 3 Group_… Group_… 0.560     2     1       0.362      0.724 0.0681  -126
## 4 Group_… Group_… 0.444     1     2       0.724      0.362 0.0807   205
## 5 Group_… Group_… 0.632     5     3       1.09       1.81  0.0590   131
## 6 Group_… Group_… 0.599     1     1       0.362      0.362 0.0839  -514
## # … with 20 more variables: lag_mm <dbl>, length <dbl>, length_mm <dbl>,
## #   matches <dbl>, matches_per_mm <dbl>, mismatches <dbl>,
## #   mismatches_per_mm <dbl>, non_cms <dbl>, non_cms_per_mm <dbl>,
## #   overlap <dbl>, rough_cor <dbl>, sd_D <dbl>, sum_peaks <dbl>,
## #   Set <chr>, BarrelA <chr>, BarrelB <chr>, BulletA <chr>, BulletB <chr>,
## #   LandA <chr>, LandB <chr>
comparisons_1 %>% 
  filter(BarrelA == "KA" & BarrelB == "KA") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("ccf scores between same barrels 1")

comparisons_1 %>% 
  filter(BarrelA == "KA" & BarrelB == "KC") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("ccf scores between different barrels 1")

comparisons_2 %>% 
  filter(BarrelA == "KC" & BarrelB == "KC") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
      facet_grid(BulletB~BulletA, labeller = "label_both") +
        xlab("Land A") +
        ylab("Land B") +
          theme(aspect.ratio = 1)+
            ggtitle("ccf scores between same barrels 2")

comparisons_2 %>% 
  filter(BarrelA == "KC" & BarrelB == "KD") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("ccf scores between different barrels 2")

comparisons_3 %>% 
  filter(BarrelA == "KF" & BarrelB == "KF") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("ccf scores between same barrels 3")

comparisons_3 %>% 
  filter(BarrelA == "KI" & BarrelB == "KJ") %>%
  ggplot(aes(x = LandA, y = LandB, fill = ccf)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = 0.5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("ccf scores between different barrels 3")

# Obtaining Random Forest Scores 
comparisons_1$rfscore <- predict(rtrees, newdata = comparisons_1, type = "prob")[,2]
comparisons_2$rfscore <- predict(rtrees, newdata = comparisons_2, type = "prob")[,2]
comparisons_3$rfscore <- predict(rtrees, newdata = comparisons_3, type = "prob")[,2]
comparisons_1 %>%
  filter(BarrelA == "KA" & BarrelB == "KA") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between same barrels 1")

comparisons_1 %>%
  filter(BarrelA == "KA" & BarrelB == "KC") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between different barrels 1")

comparisons_2 %>%
  filter(BarrelA == "KC" & BarrelB == "KC") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between same barrels 2")

comparisons_2 %>%
  filter(BarrelA == "KC" & BarrelB == "KD") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between different barrels 2")

comparisons_3 %>%
  filter(BarrelA == "KI" & BarrelB == "KI") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between same barrels 3")

comparisons_3 %>%
  filter(BarrelA == "KI" & BarrelB == "KJ") %>%
  ggplot(aes(x = LandA, y = LandB, fill = rfscore)) +
    geom_tile() +
      scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
        facet_grid(BulletB~BulletA, labeller = "label_both") +
          xlab("Land A") +
          ylab("Land B") +
            theme(aspect.ratio = 1)+
              ggtitle("rf scores between different barrels 3")

Bullet Scores

Currently on hold

Bullet_Scores_1 <- comparisons_1 %>% 
  group_by(BulletA, BulletB) %>% 
    tidyr::nest()

Bullet_Scores_1 <- Bullet_Scores_1 %>% 
  mutate(Bullet_Score = data %>% 
          future_map_dbl(.f = function(d) max(compute_average_scores(land1 = d$LandA, land2 = d$LandB, d$rfscore))))

Bullet_Scores_1 %>% 
  select(-data) %>% 
    arrange(desc(Bullet_Score))
## # A tibble: 121 x 3
##    BulletA BulletB Bullet_Score
##    <chr>   <chr>          <dbl>
##  1 U10     U10            0.977
##  2 U15     U15            0.977
##  3 U28     U28            0.977
##  4 U36     U36            0.977
##  5 U37     U37            0.977
##  6 U40     U40            0.977
##  7 U42     U42            0.977
##  8 U77     U77            0.977
##  9 U40     U36            0.708
## 10 U36     U40            0.708
## # … with 111 more rows
Bullet_Scores_1 <- Bullet_Scores_1 %>% 
  mutate(data = data %>% 
           future_map(.f = function(d){
              d$samepath = bullet_to_land_predict(land1 = d$LandA, land2 = d$LandB, d$rfscore, difference = 0.1)
  d
}))

Bullet_Scores_1 %>% 
  tidyr::unnest(data) 
## # A tibble: 19,044 x 32
##    BulletA BulletB Bullet_Score Bullet1 Bullet2   ccf   cms  cms2
##    <chr>   <chr>          <dbl> <chr>   <chr>   <dbl> <dbl> <dbl>
##  1 1       1              0.464 Group_… Group_… 1        27    14
##  2 1       1              0.464 Group_… Group_… 0.566     1     2
##  3 1       1              0.464 Group_… Group_… 0.474     3     3
##  4 1       1              0.464 Group_… Group_… 0.552     2     2
##  5 1       1              0.464 Group_… Group_… 0.400     0     0
##  6 1       1              0.464 Group_… Group_… 0.385     5     4
##  7 1       1              0.464 Group_… Group_… 0.381     3     3
##  8 1       1              0.464 Group_… Group_… 0.554     1     1
##  9 1       1              0.464 Group_… Group_… 0.477     1     1
## 10 1       1              0.464 Group_… Group_… 0.303     2     1
## # … with 19,034 more rows, and 24 more variables: cms2_per_mm <dbl>,
## #   cms_per_mm <dbl>, D <dbl>, lag <dbl>, lag_mm <dbl>, length <dbl>,
## #   length_mm <dbl>, matches <dbl>, matches_per_mm <dbl>,
## #   mismatches <dbl>, mismatches_per_mm <dbl>, non_cms <dbl>,
## #   non_cms_per_mm <dbl>, overlap <dbl>, rough_cor <dbl>, sd_D <dbl>,
## #   sum_peaks <dbl>, Set <chr>, BarrelA <chr>, BarrelB <chr>, LandA <chr>,
## #   LandB <chr>, rfscore <dbl>, samepath <lgl>
Bullet_Scores_2 <- comparisons_2 %>% 
  group_by(BulletA, BulletB) %>% 
    tidyr::nest()

Bullet_Scores_2 <- Bullet_Scores_2 %>% 
  mutate(Bullet_Score = data %>% 
          future_map_dbl(.f = function(d) max(compute_average_scores(land1 = d$LandA, land2 = d$LandB, d$rfscore))))

Bullet_Scores_2 %>% 
  select(-data) %>% 
    arrange(desc(Bullet_Score))
## # A tibble: 121 x 3
##    BulletA BulletB Bullet_Score
##    <chr>   <chr>          <dbl>
##  1 U23     U23            0.977
##  2 U34     U34            0.977
##  3 U41     U41            0.977
##  4 U61     U61            0.977
##  5 U63     U63            0.977
##  6 U66     U66            0.977
##  7 U73     U73            0.977
##  8 U98     U98            0.977
##  9 U98     U61            0.546
## 10 U61     U98            0.546
## # … with 111 more rows
Bullet_Scores_2 <- Bullet_Scores_2 %>% 
  mutate(data = data %>% 
           future_map(.f = function(d){
              d$samepath = bullet_to_land_predict(land1 = d$LandA, land2 = d$LandB, d$rfscore, difference = 0.1)
  d
}))

Bullet_Scores_2 %>% 
  tidyr::unnest(data) 
## # A tibble: 19,044 x 32
##    BulletA BulletB Bullet_Score Bullet1 Bullet2   ccf   cms  cms2
##    <chr>   <chr>          <dbl> <chr>   <chr>   <dbl> <dbl> <dbl>
##  1 1       1              0.460 Group_… Group_… 1        27    13
##  2 1       1              0.460 Group_… Group_… 0.536     1     0
##  3 1       1              0.460 Group_… Group_… 0.560     2     1
##  4 1       1              0.460 Group_… Group_… 0.444     1     2
##  5 1       1              0.460 Group_… Group_… 0.632     5     3
##  6 1       1              0.460 Group_… Group_… 0.599     1     1
##  7 1       1              0.460 Group_… Group_… 0.588     2     1
##  8 1       1              0.460 Group_… Group_… 0.400     1     1
##  9 1       1              0.460 Group_… Group_… 0.358     3     2
## 10 1       1              0.460 Group_… Group_… 0.312     3     2
## # … with 19,034 more rows, and 24 more variables: cms2_per_mm <dbl>,
## #   cms_per_mm <dbl>, D <dbl>, lag <dbl>, lag_mm <dbl>, length <dbl>,
## #   length_mm <dbl>, matches <dbl>, matches_per_mm <dbl>,
## #   mismatches <dbl>, mismatches_per_mm <dbl>, non_cms <dbl>,
## #   non_cms_per_mm <dbl>, overlap <dbl>, rough_cor <dbl>, sd_D <dbl>,
## #   sum_peaks <dbl>, Set <chr>, BarrelA <chr>, BarrelB <chr>, LandA <chr>,
## #   LandB <chr>, rfscore <dbl>, samepath <lgl>
Bullet_Scores_3 <- comparisons_3 %>% 
  group_by(BulletA, BulletB) %>% 
    tidyr::nest()

Bullet_Scores_3 <- Bullet_Scores_3 %>% 
  mutate(Bullet_Score = data %>% 
          future_map_dbl(.f = function(d) max(compute_average_scores(land1 = d$LandA, land2 = d$LandB, d$rfscore))))

Bullet_Scores_3 %>% 
  select(-data) %>% 
    arrange(desc(Bullet_Score))
## # A tibble: 121 x 3
##    BulletA BulletB Bullet_Score
##    <chr>   <chr>          <dbl>
##  1 U14     U14            0.977
##  2 U27     U27            0.977
##  3 U33     U33            0.977
##  4 U36     U36            0.977
##  5 U45     U45            0.977
##  6 U49     U49            0.977
##  7 U56     U56            0.977
##  8 U65     U65            0.977
##  9 U36     U33            0.782
## 10 U33     U36            0.782
## # … with 111 more rows
Bullet_Scores_3 <- Bullet_Scores_3 %>% 
  mutate(data = data %>% 
           future_map(.f = function(d){
              d$samepath = bullet_to_land_predict(land1 = d$LandA, land2 = d$LandB, d$rfscore, difference = 0.1)
  d
}))

Bullet_Scores_3 %>% 
  tidyr::unnest(data) 
## # A tibble: 19,044 x 33
##    BulletA BulletB Bullet_Score Bullet1 Bullet2 Legacy_Features   ccf   cms
##    <chr>   <chr>          <dbl> <chr>   <chr>   <list>          <dbl> <dbl>
##  1 1       1              0.444 Group_… Group_… <df[,17] [1 × … 1        34
##  2 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.301     1
##  3 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.363     2
##  4 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.330     2
##  5 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.361     3
##  6 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.370     3
##  7 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.518     1
##  8 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.362     1
##  9 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.320     2
## 10 1       1              0.444 Group_… Group_… <df[,17] [1 × … 0.155     3
## # … with 19,034 more rows, and 25 more variables: cms2 <dbl>,
## #   cms2_per_mm <dbl>, cms_per_mm <dbl>, D <dbl>, lag <dbl>, lag_mm <dbl>,
## #   length <dbl>, length_mm <dbl>, matches <dbl>, matches_per_mm <dbl>,
## #   mismatches <dbl>, mismatches_per_mm <dbl>, non_cms <dbl>,
## #   non_cms_per_mm <dbl>, overlap <dbl>, rough_cor <dbl>, sd_D <dbl>,
## #   sum_peaks <dbl>, Set <chr>, BarrelA <chr>, BarrelB <chr>, LandA <chr>,
## #   LandB <chr>, rfscore <dbl>, samepath <lgl>
Bullet_Scores_Examin_1 <- Bullet_Scores_1 %>% 
  tidyr::unnest(data)

#Bullet_Scores_Examin_1


Bullet_Scores_Examin_2 <- Bullet_Scores_2 %>% 
  tidyr::unnest(data)

#Bullet_Scores_Examin_2


Bullet_Scores_Examin_3 <- Bullet_Scores_3 %>% 
  tidyr::unnest(data)

#Bullet_Scores_Examin_3
ggplot(subset(Bullet_Scores_Examin_1, !BulletB %in% c("1", "2", "3")))+
      geom_tile(aes(x = BulletA, y = BulletB, fill = rfscore))+
        scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
  facet_grid(BarrelB~BarrelA, labeller = "label_both", scales = "free")+
    ggtitle("Group1")

ggplot(subset(Bullet_Scores_Examin_2, !BulletB %in% c("1", "2", "3")))+
      geom_tile(aes(x = BulletA, y = BulletB, fill = rfscore))+
        scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
  facet_grid(BarrelB~BarrelA, labeller = "label_both", scales = "free")+
    ggtitle("Group2")

ggplot(subset(Bullet_Scores_Examin_3, !BulletB %in% c("1", "2", "3")))+
      geom_tile(aes(x = BulletA, y = BulletB, fill = rfscore))+
        scale_fill_gradient2(low = "grey80", high = "darkorange", 
                       midpoint = .5) +
  facet_grid(BarrelB~BarrelA, labeller = "label_both", scales = "free")+
    ggtitle("Group3")